Nature-inspired microfluidic propulsion using magnetic actuation 
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In this work we mimic the efficient propulsion mechanism of natural cilia by magnetically actuating 
thin films in a cyclic but non-reciprocating manner. By simultaneously solving the elasto-dynamic, 
magnetostatic and fluid mechanics equations, we show that the amount of fluid propelled is propor- 
tional to the area swept by the cilia. By using the intricate interplay between film magnetization 
and applied field we are able to generate a pronounced asymmetry and associated flow. We delineate 
the functional response of the system in terms of three dimensionless parameters that capture the 
relative contribution of elastic, inertial, viscous and magnetic forces. 
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A rapidly growing field in biotechnology is the use of 
lab-on-a-chip devices to analyse bio-fluids [J, 0, ■ Such 
fluids have to be preprocessed (for example, mixed with 
other fluids jj]) and transported to and from one or many 
micro-chambers where the biochemical analyses are per- 
formed. The microfiuid transport through these stages is 
usually performed by downscaling conventional methods 
such as syringe pumps [H, @] , micropumps Q , or by ex- 
ploiting electro-magnetic actuation, as in electro-osmotic 
0, [1] and magnetohydrodynamic devices 0, [13]. How- 
ever, when transporting biological fluids (which usually 
have high conductivity), the use of electric fields may in- 
duce heating, bubble formation and pH gradients from 

In this work, we 
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explore a new way to manipulate fluids in microfluidic 
systems, inspired by nature, through the magnetic actu- 
ation of artifical cilia. 

Fluid dynamics at the micrometer scale is dominated 
by viscosity rather than inertia. This has important 
consequences for fluid propulsion mechanisms [141 ] - In 
particular, mechanical actuation will only be effective in 
propelling fluids if their motion is cyclic, but asymmet- 
ric in shape change. Nature has solved this problem by 
means of hair-like structures, called cilia, whose beating 
pattern is asymmetric and consists of an effective and a 
recovery stroke [HJ. While natural cilia use an internal 
forcing system based on motor proteins (dyneins), the 
key challenge for its artificial equivalent is the design of 
an externally-applied loading system that will generate a 
similar non-reciprocating motion. Recently, electrostatic 
artificial cilia have been experimentally shown to induce 
effective micro- mixing Q. In addition, magnetic fields 
are also used to induce flow, but the asymmetry gener- 
ated was found to be relatively small [l6| . In this work we 
report on the identification of two simple magnetically- 
driven configurations that can create a large asymme- 
try. We will show that the fluid propelled is linearly 
proportional to the swept area by the film (the configu- 



rational space), which has been shown so far only for a 
non-actuated kinetic three-sphere model 17|] . The first is 
based on a magnetic instability that develops when the 
applied magnetic field is opposite to the direction of the 
magnetization in a permanently magnetic film. In a sec- 
ond configuration we will demonstrate that asymmetry 
can be achieved in a super-paramagnetic film, based on 
the intricate inter-play between the geometry of the film, 
the externally-applied field and the internally-induced 
magnetization. 

The numerical model used in this work is based on 
a two-dimensional finite clement representation of thin 
magnetic films, employing Euler-Bcrnoulli beam ele- 
ments. We simultaneously solve for the elasto-dynamic 
equations of motion and Maxwell's equations, so that we 
can accurately account for the elastic, inertia and mag- 
netic interaction in a non-linear geometry setting. We 
explicitly couple this Lagrangian solid-dynamics model 
to an Eulcrian fluid dynamics model through Lagrange 
multipliers. Input to the fluid mechanics model are the 
positions and velocities of the film at all times which re- 
sult in a full velocity field in the fluid. We calculate the 
drag forces on the film as tractions via the stress tensor 
in the fluid. The traction distribution is subsequently 
imposed as surface tractions in the magneto-mechanical 
model. A detailed description of the model can be found 
in the Supporting Information [l8j ]. 

We study a periodic arrangement of permanently mag- 
netic (PM) cilia in a microfluidic channel of height 5L, 
with the cilia spaced 5L apart, where L is the length 
of the cilia. A square unit-cell is identified consisting 
of one cilium. No-slip boundary conditions are applied 
at the top and bottom boundaries of the channel and 
periodic boundary conditions at the left and right ends 
of the unit-cell. The fluid has a viscosity /i = 1 mPas. 
The film has a thickness h = 2 /jm, effective stiffness 
E = E/(l — v 2 ) ~ 1 MPa, where E is the elastic modulus 
and v is the Poisson's ratio, and density p = 1600 kg/m 3 . 
The initial geometry of the film is a quarter of a circle 
with radius 100 fiia fixed at the bottom of the channel, 
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see instant 1 in Fig. 1(a) The direction of the magne- 



2 




(a) 



0.2 0.4 0.6 0.8 1 
Normalised coordinate along the film (Q 



(b) 




Fixed end 



(a) 



0.2 0.4 0.6 0.8 1 
Normalised coordinate along the film (Q 



(b) 



FIG. 1: Buckling of a curled permanently magnetic (PM) 
film as a result of magnetic actuation, during the propulsion 
of fluid, (a) Snapshots of the film at ms, 0.3 ms, 0.6 ms, 1.1 
ms and 3 ms. (b) Normalized torque distribution along the 
film corresponding to the snapshots shown in (a). 



tization is along the film with the magnetization vector 
pointing from the fixed end to the free end. The remnant 
magnetization of the film is taken to be M r = 15 kA/m. 
A uniform external field of magnitude Bq = 13.3 mT is 
applied at 225° to the x axis from t = ms to t = 1 ms 
and then linearly reduced to zero in the next 0.2 ms. The 
results of the non-reciprocating motion of the film in the 



fluid during magnetic actuation are shown in Fig. 1(a) 



The Eulerian fluid mesh is not shown for clarity. When 
the external field is applied, clockwise torques (N z is the 
magnetic body torque) are acting on the portion near 
the fixed end of the film while near the free end counter- 



clockwise torques develop (see instance 1 in Fig. 1(b) I 



Under the influence of such a system of moments, the 
film undergoes a buckling kind of instability. This can 



be nicely seen from instances 1 and 2 in Figs. 1(a) and 



1(b) During this stage the position of zero torque is al- 



most fixed, while the torques at the free end increase. 
This causes the film to snap through to configurations 
3 and 4 during which the zero-torque position travels to 
the fixed end. Clearly, the initially opposing directions 
of the internal magnetization and the applied magnetic 
field are essential in generating an instability that causes 
a large bending deformation during application of the 
field. Then, the applied field is reduced to zero and the 
film returns to the initial position through instance 5 in 



Fig. 1 (a) Note that the propulsive action in the effective 



stroke (red) takes place during the elastic recovery of the 
film, while the film stays low in the recovery stroke due 
to the buckling-enforced snap-through. 

For a PM film the torques are maximum when the lo- 
cal magnetic film is perpendicular to the (remnant) mag- 
netization. For a super-paramagnetic (SPM) film, how- 
ever, the magnetization is induced by the field itself, pos- 
ing different requirements on the applied magnetic fields 
in order to deform the film. A straight, magnetically 
anisotropic SPM film (having susceptibilities 4.6 and 0.8 
in the tangential and normal directions, respectively), is 
subjected to a magnetic field with magnitude Bq = 31.5 



FIG. 2: Motion of a super-paramagnetic (SPM) film in a 
rotating magnetic field, during the propulsion of fluid, (a) 
Snapshots of the film at ms, 2.5 ms, 5.0 ms, 7.5 ms and 
8.5 ms. (b) Normalized torque distribution along the film 
corresponding to the snapshots shown in (a). Here fio is the 
permeability of vacuum, h is the thickness at position £ and 
h( = o is the thickness at the fixed end. 



mT that is rotated from 0° to 180° in t = 10 ms and 
then kept constant during the rest of the cycle. The film 
has a length L = 100 fim, effective stiffness E = 1 MPa 
and density p = 1600 kg/m 3 . Its cross-section is tapered, 
with the thickness varying linearly along its length, hav- 
ing h = 2 [im at the left (attached) end and h = 1 fj,m 
at the right end. Figure 2(a) shows that in the effective 
stroke the portion of the beam near the free end is nearly 
straight. This is due to the fact that in this region the 
film can easily follow the applied field so that field and 
magnetization arc almost parallel, causing the magnetic 
torque to be low in this region of the film (instances 2, 
3 and 4 in Fig. 



2(b)). When the film has reached po- 



sition 4, the magnetization in the film is such that the 
torques are oriented clockwise near the fixed end and an- 
ticlockwise near the free end, resulting in strong bending 
of the film. From Fig. |2(b)| it can be seen that during 
the recovery stroke (in black) the position of zero torque 
propagates from the fixed end to the free end (from in- 
stance 4 to 5). Here the tapering is essential, causing the 
torque per unit length to be higher at the fixed end, al- 
lowing the film to recover to the initial position (1). This 
behaviour is very similar to that of natural cilia [151 ]. It 
is to be noted that the film recovers in the presence of an 
applied magnetic field. This sensitive interplay between 
stored elastic energy and controlled applied field can be 
exploited to provide a large asymmetry in motion. 

Next we analyze how much fluid is propelled by the two 
cases analyzed. We record the fluid volume transported 
through the channel per cycle and per unit out-of-planc 
thickness, giving an area flow per cycle. As a measure 
for the asymmetry, we compute the area swept by the 
free end of the film during one cycle (i.e. the area en- 
closed by the dashed lines in Figs. 1(a) and 2(a)) and 
vary this area by tuning the magnitude of the applied 
magnetic field (all other parameters remain unchanged). 
Fig. [3] shows the area flow per cycle as a function of the 
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swept area for several different cases. We have normal- 
ized both quantities by the maximum area that the tip 
can sweep, ttL 2 /2. For three values of the magnetic field 
we plot the film tip trajectories for the PM and SPM 
configurations. The cycle times are 35 and 10 ms, re- 
spectively. The flux across the channel shows a linear 
dependence on the swept area. Similar result has been 
shown in where it is shown that the velocity of a 
three sphere swimmer is proportional to the area swept 
in the configurational space. Due to the linear correla- 
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Normalized swept area 
FIG. 3: Variation of normalized area flow with swept area. 

tion between the swept area and the fluid flow, the swept 
area can be used as a measure of effectiveness of the ac- 
tuator, representing the fluid volume displaced. This al- 
lows uncoupling the magneto-mechanical motion of the 
cilia from the computationally-intensive fluid dynamics 
calculations. Instead, we account for the fluid by means 
of velocity-proportional drag forces (using resistive force 
theory [20j) on the cilia, with the drag coefficients cali- 
brated to the coupled solid-fluid model (see (HI). 

To identify the dimensionlcss parameters that govern 
the behavior of the system, we start from the virtual 
work equation for the film [18[, neglecting the axial de- 
formations: J EIv"Sv"dx + f pAvSvdx — J N z Sv'Adx + 
J CyiSvbdx = 0, where ( )' = d( )/dx, (') = d{ )/dt, 
I = bh 3 I '12 is the second moment of area with b the 
out-of-plane thickness, A = bh is the cross-sectional area 
of the film and v is the transverse displacement. In the 
virtual work equation the first, second, third and last 
terms respectively represent the virtual work done by the 
clastic internal bending moments, the inertial forces, the 
magnetic couple and the fluid drag forces. We intro- 
duce the dimensionlcss variables V, T and X, such that 
v = VL, x = XL and t = Tt rc {, where L is a characteris- 
tic length (taken to be the length of the film) and i re f a 
characteristic time. Substitution of these variables in the 
virtual work equation and normalization with the elas- 
tic term reveals the three governing dimensionless num- 
bers: the inertia number, /„ = 12pL 4 / Eh 2 t 2 . cf , i.e. the 
ratio of inertial to elastic force, the magnetic number 
M n = 12N Z L 2 1 Eh 2 , i.e. the ratio of magnetic to elastic 
force and the fluid number F n = l2C y L 4 / Eh 3 t rc i, the 
ratio of fluid to elastic force. By substituting the torque 



expression for the two different magnetic materials, the 
magnetic number M n for the PM film is linear in the ap- 
plied field, l2M r B L 2 /Eh 2 , while for the SPM film it is 
quadratic, 12BqL 2 / fi Eh 2 . 

We proceed by exploring the functional response of 
the system in terms of the swept area and cycle time, in 
dependence of the three dimensionless parameters. We 
analyzed many different combinations of /„ , M n and F n , 
the results of which are summarized in Figs. 4(a)|4(c)" 
for the PM sy stem and in Fi gs. |4(d)||4(f)| for the SPM 
system. Figs. 4(a) and |4(d)| show the swept area as a 
function of M n for several combinations of /„ and F n . 
The combinations are indicated by the different symbols, 
corresponding to specific locations in Figs. |4(b)| and 4(e) 



The effect of all three parameters can also be nicely sum- 
marised by analyzing what magnetic number and cycle 
time is needed to sweep a normalized area of 0.2, for a 
given range of /„ and F n values (see Figs. |4(b)| |4(c)| 
4(c)| |4(f)[ ) . The swept area increases with M n reaching a 



maximum of 0.4 for the PM system (see Fig. |4(a)[ ), while 
values of 0.7 can be reached by the SPM system (see 
Fig. |4(d)[ ). For both systems the M n needed strongly 
increases with F n . In other words, for a given elastic 
parameter set, larger magnetic forces are needed to over- 



4(b) 



4(d) 



come the drag forces imposed by the fluid (see F igs. |4(a) 
|4(d)| and 4(e) I. It can be seen from Figs. 4(a) and 
that the effect of F n is gradual for the PM sys 



tern, while for the SPM case it is absent for small F n , 
but suddenly kicks in for F n larger than 10. In addi- 
tion, the inertial forces assist in generating asymmetry 
for both cases, although for the SPM system inertial ef- 
fects are only triggered for very large fluid numbers (see 
Fig. [4(d)]). " 

For the analysis of the normalised cycle time the ref- 
erence time was taken to be the time during which the 
field was applied (PM) or rotated (SPM). It was observed 
that the cycle time dependence on M n and /„ was very 
weak and mostly completely absent. The only clear de- 
pendence found was on F n which we show in the inset 
of Figs. 4(a) and |4(d)| in terms of the normalized cy- 
cle times averaged over all different M n values analyzed 
and in Figs. 4(c)| and |4(f)| as the time required to sweep 
a normalised area of 0.2. For the PM case the effective 
stroke is generated through the clastic recovery of the de- 
formed film, without noticeable effect of inertial forces. 
For such ovcrdamped systems the time taken by the sys- 
tem to return to the initial position scales linearly with 
F n , the ratio of fluid to elastic forces. The variation of 
normalized cycle time with F n for the SPM case is much 
smaller. This is due to the fact that the total cycle is 
performed in the presence of magnetic forces. For small 
F n and /„ the mean normalized cycle time (see Fig. |4(f)[ ) 
is approximately equal to one; only for large F n and /„ 
the cycle time is increased. At large F n the system re- 
lies on the recovery (going from instance 5 to instance 1 
in Fig. \2§ of the curved tip against high viscous forces 
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FIG. 4: Functional response of the PM system (top row, [(a)f(c)j ) and the SPM system (bottom row, |(d)|(f)) . |(a)|(d)| Normalized 



swept area as a function of M n for several combinations of I n and F n , corresponding to the symbols of Figs, (b) and (e) The 
inset shows the mean normalized cycle time as a function of F„. The mean is obtained by averaging all times corresponding to 
the data points that make up the specific M n -swept area curve. |(b)|(e)| Contours of M n needed to sweep a normalized area of 
0.2 for a wide range of I„ and F n values. (c)|(f) Contours of normalized cycle time corresponding to (b) and (e) 



(see Figs. |4(d)| |4(f)[ ). The systems demonstrate an un- 
dcrdamped behaviour at large /„ values causing inertial 
forces to generate large oscillations leading to a larger 
normalised time (see Fig. |4(f)[ ). For a normalized area 
of 0.2 the response of both systems in the range /„ < 0.1 
and F n < 10 is quasi-static, i.e. independent of inertial 
and viscous effects. 

To summarize, we have proposed and analysed mag- 
netic artificial cilia which can transport fluid in microflu- 
idic channels. The main result is that we have found two 
simple and novel actuation mechanisms which can gen- 
erate a pronounced asymmetric motion of the cilia. One 
configuration is based on the buckling of a permanently 
magnetic film and the other is based on the intricate in- 
teraction between the applied field and the magnetization 
in a super-paramagnetic film. We have shown that the 
fluid propelled is linearly proportional to the area swept 
by the film, which has so far only been shown for a non- 
actuated kinetic system [17]. Finally, we have identified 



the range of dimcnsionlcss parameters for which the ar- 
tificial cilia exhibit an optimal behavior. The analysis 
presented can be used as a guideline to make artificial 
cilia for microfluidic transport in lab-on-a-chip systems. 

This work is a part of the 6 th Framework European 
project 'Artie', under contract STRP 033274. 
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